### Install Packages:
rm(list = ls())
library(tidyverse)
library(lmtest)
library(sandwich)
library(ggplot2)
library(ggthemes)
library(ggpubr)
library(stargazer)
library(cregg)
library(dotwhisker)
library(dplyr)
### Set Working Directory
# setwd()
source("functions/vcovCluster.R")
source("functions/effect_plotter.R")
### Load Data
conjoint = readRDS("religion_conjoint_reshaped.rds")
head(conjoint)
### Cleaning Attributes
conjoint = conjoint %>%
mutate(trust = recode(trust, "Not trust at all" = "0", "A little" = "1", "Somewhat" = "2", "A lot" = "3",
"Completely trust" = "4", .default = NA_character_),
trust = as.numeric(trust),
cjt_nationality = recode(cjt_nationality, "$q_country" = "Home"),
cjt_age = factor(cjt_age, levels = c("28", "45", "70")),
cjt_education = factor(cjt_education, levels = c("Self-taught in religion", "BA in Engineering", "BA in Islamic Fiqh",
"Doctorate in Engineering", "Doctorate in Islamic Fiqh")),
cjt_nationality = factor(cjt_nationality, levels = c("Home", "Iraq", "Saudi Arabia")),
cjt_hafiz = factor(cjt_hafiz, levels = c("No", "Yes")),
cjt_location = factor(cjt_location, levels = c("Internet", "Government Television", "Satellite Television",
"Mosque registered with government", "Mosque not registered with government")),
cjt_followers = factor(cjt_followers, levels = c("5,000", "50,000", "500,000")),
cjt_ideology = factor(cjt_ideology, levels = c("Sunni Orthodox", "Muslim Brotherhood", "Salafi",
"Sufi", "Shiite")),
cjt_usa = factor(cjt_usa, levels = c("No Opinion", "Participates in meetings with US embassy",
"Criticizes United States", "Supports armed resistance against United States")),
cjt_politics = factor(cjt_politics, levels = c("None", "Supports Government", "Criticizes Government")))
age = unique(conjoint$cjt_age)
education = unique(conjoint$cjt_education)
nationality = unique(conjoint$cjt_nationality)
hafiz = unique(conjoint$cjt_hafiz)
location = unique(conjoint$cjt_hafiz)
followers = unique(conjoint$cjt_followers)
ideology = unique(conjoint$cjt_ideology)
usa = unique(conjoint$cjt_usa)
politics = unique(conjoint$cjt_politics)
refcats = data_frame(variable = c("28", "Self-taught in religion", "Home", "No", "Internet", "5,000",
"Sunni Orthodox", "No Opinion", "None"),
coef = 0, se = 0, attribute = 0)
attributes = data_frame(variable = c("Age", "Education", "Hafiz", "Nationality", "Location",
"Followers", "Ideology", "USA", "Politics"),
coef = NA_real_,
se = NA_real_, attribute = 1)
# Subset data to Sunni Muslims
conjoint = conjoint %>% filter(sect == "Sunni"|sect=="Salafi")
### Primary Conjoint Results: Figure 1 from Article
df = conjoint
# Regression Model
basemod1 = lm(binary_trust ~  cjt_age + cjt_education + cjt_hafiz + cjt_nationality + cjt_location +
cjt_followers + cjt_ideology + cjt_usa + cjt_politics, data = df)
mod1<-summary(basemod1)
basemod1 = coeftest(basemod1, vcov = vcovCluster(basemod1, factor(df$id)))
basemod1
# Robustness Check with Fixed Effects for Respondent Country
fe_mod = lm(binary_trust ~  cjt_age + cjt_education + cjt_hafiz + cjt_nationality + cjt_location +
cjt_followers + cjt_ideology + cjt_usa + cjt_politics + factor(residence_country), data = df)
fe_mod = coeftest(fe_mod, vcov = vcovCluster(fe_mod, factor(df$id)))
fe_mod
stargazer(fe_mod)
# Subsetting to Islamists
df = conjoint %>% filter(islamist==1)
# Primary Results Among Islamists (coefficients reported in Figure 2 of Article)
basemod = lm(binary_trust ~  cjt_age + cjt_education + cjt_hafiz + cjt_nationality + cjt_location +
cjt_followers + cjt_ideology + cjt_usa + cjt_politics, data = df)
basemod = coeftest(basemod, vcov = vcovCluster(basemod, factor(df$id)))
basemod
plot_islamist = effect_plotter(model.output =  basemod, names.variables = c("Age", "Education", "Hafiz", "Nationality", "Location",
"Followers", "Ideology", "USA", "Politics"),
names.levels = list(c("28", "45", "70"), c("Self-taught in religion", "BA in Engineering",
"BA in Islamic Fiqh","Doctorate in Engineering",
"Doctorate in Islamic Fiqh"),c("No", "Yes"),
c("Home", "Iraq", "Saudi Arabia"),
c("Internet", "Government Television", "Satellite Television",
"Mosque registered w/ gov", "Mosque not registered w/ gov"),
c("5,000", "50,000", "500,000"), c("Sunni Orthodox", "Muslim Brotherhood", "Salafi",
"Sufi", "Shiite"),
c("No opinion", "Meets with US embassy",
"Criticizes US", "Supports armed resistance against US"),
c("None", "Supports government", "Criticizes government")),
effect.label = "Change in Probability",
x.lower = -0.4, x.upper = 0.2, labs = T, title = "Islamists")
# Subsetting to Non-Islamists
df = conjoint %>% filter(islamist==0)
# Primary Results Among Non-Islamists (coefficients reported in Figure 2 of Article)
basemod = lm(binary_trust ~  cjt_age + cjt_education + cjt_hafiz + cjt_nationality + cjt_location +
cjt_followers + cjt_ideology + cjt_usa + cjt_politics, data = df)
basemod = coeftest(basemod, vcov = vcovCluster(basemod, factor(df$id)))
basemod
plot_not_islamist = effect_plotter(model.output =  basemod, names.variables = c("Age", "Education", "Hafiz", "Nationality", "Location",
"Followers", "Ideology", "USA", "Politics"),
names.levels = list(c("28", "45", "70"), c("Self-taught in religion", "BA in Engineering",
"BA in Islamic Fiqh","Doctorate in Engineering",
"Doctorate in Islamic Fiqh"),c("No", "Yes"),
c("Home", "Iraq", "Saudi Arabia"),
c("Internet", "Government Television", "Satellite Television",
"Mosque registered w/ gov", "Mosque not registered w/ gov"),
c("5,000", "50,000", "500,000"), c("Sunni Orthodox", "Muslim Brotherhood", "Salafi",
"Sufi", "Shiite"),
c("No opinion", "Meets with US embassy",
"Criticizes US", "Supports armed resistance against US"),
c("None", "Supports government", "Criticizes government")),
effect.label = "Change in Probability",
x.lower = -0.4, x.upper = 0.2, labs = F, title = "Non-Islamists")
#ggarrange(p1, p2, widths = c(3.55, 1.45))%>%
islamist_plots<-ggarrange(plot_islamist,plot_not_islamist,ncol=2,widths = c(3.05, 1.6))
islamist_plots
term<-c("Sunni Orthodox","Muslim Brotherhood","Salafi","Sufi","Shiite","Sunni Orthodox","Muslim Brotherhood","Salafi","Sufi","Shiite")
estimate<-c(0,-0.094,-0.109,-0.203,-0.396,0,-0.161,-0.124,-0.140,-0.343)
std.error<-c(0,0.014,0.014,0.014,0.013,0,0.009,0.009,0.009,0.009)
model<-c("Islamist","Islamist","Islamist","Islamist","Islamist","Non-Islamist","Non-Islamist","Non-Islamist","Non-Islamist","Non-Islamist")
islamist_df<-cbind.data.frame(term,estimate,std.error,model)
islamist_df
islamist_df$term<-factor(islamist_df$term,levels=c("Sunni Orthodox","Muslim Brotherhood","Salafi","Sufi","Shiite"))
islamist_subgroup_plot<-dwplot(islamist_df, dot_args = list(size = 2.5),whisker_args = list(size = 1)) +
theme_bw() + #geom_point(size = 3, shape = 16) +
geom_vline(xintercept = 0, colour = "black", linetype = 2) +
ggtitle("") +
scale_x_continuous("AMCE",limits = c(-.45, .05)) + theme(panel.spacing = unit(1, "lines")) +
scale_y_discrete("Ideology Levels\n",labels=c("Shiite","Sufi","Salafi","Muslim \nBrotherhood","Sunni \nOrthodox")) +
theme(axis.title=element_text(size=10),
legend.text=element_text(size=10),legend.position="bottom",
axis.text.x=element_text(size=10,color="black"),
axis.text.y=element_text(size=10,color="black")) + scale_color_manual(values=c("darkgrey", "black")) +
labs(color='Respondent Ideology')
islamist_subgroup_plot
